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Abstract — In this paper, we introduce an analytical framework 
to compute the average rate of downlink heterogeneous cellular 
networks. The framework leverages recent application of stochas- 
tic geometry to other-cell interference modeling and analysis. The 
heterogeneous cellular network is modeled as the superposition of 
many tiers of Base Stations (BSs) having different transmit power, 
density, path-loss exponent, fading parameters and distribution, 
and unequal biasing for flexible tier association. A long-term 
averaged maximum biased-received-power tier association is 
considered. The positions of the BSs in each tier are modeled 
as points of an independent Poisson Point Process (PPP). Under 
these assumptions, we introduce a new analytical methodology 
to evaluate the average rate, which avoids the computation of 
the Coverage Probability (Pcov) and needs only the Moment 
Generating Function (MGF) of the aggregate interference at 
the probe mobile terminal. The distinguishable characteristic of 
our analytical methodology consists in providing a tractable and 
numerically efficient framework that is applicable to general 
fading distributions, including composite fading channels with 
small- and mid-scale fluctuations. In addition, our method can 
efflciently handle correlated Log-Normal shadowing with little 
increase of the computational complexity. The proposed MGF- 
based approach needs the computation of either a single or 
a two-fold numerical integral, thus reducing the complexity 
of Pcov-based frameworks, which require, for general fading 
distributions, the computation of a four-fold integral. 

Index Terms — Heterogeneous cellular networks, aggregate in- 
terference modeling, stochastic geometry, average rate. 



I. Introduction 

THE analytical performance modeling of cellular networks 
is a long-standing open issue [1]. The mathematical 
intractability mainly arises from the difficulty of accurately 
modeling other-cell interference by taking into account the 
spatial positions of the Base Stations (BSs) and the stochastic 
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character of the wireless channel [2]. For this reason, accurate 
performance analysis is usually conducted via costly, time- 
consuming, and often proprietary system-level simulators [3], 
This approach, however, seldom provides insightful informa- 
tion on system design and on the dependency of the system 
parameters to optimize. This situation is even exacerbated 
in future cellular deployments, which are becoming more 
heterogeneous with the introduction of new infrastructure 
elements, e.g., femto/pico BSs, fixed/mobile relays, cognitive 
radios, and distributed antennas [4]-[9] and [10] for a survey. 

A. Abstraction Models for Analysis and Design of Cellular 
Networks 

To circumvent this problem, communications theorists usu- 
ally resort to "abstractions" for tractable other-cell interference 
modeling and for performance analysis. These abstractions 
usually encompass simplified spatial models for the locations 
of the BSs. In particular, three abstraction models are com- 
monly adopted: i) the Wyner model [11]; ii) the single-cell 
interfering model [12]; and iii) the regular hexagonal or square 
grid model [13]. These abstraction models, however, are often 
either over-simplistic or inaccurate [14]. Furthermore, in some 
cases, as for the regular hexagonal/square grid model, they 
still require either intensive numerical simulations or multi- 
fold numerical integrations. Finally, these abstraction models 
usually provide information for specific BSs deployments, and 
typically fail to provide useful information for more random, 
unplanned, and/or clustered BSs deployments, which are typ- 
ical of emerging heterogeneous cellular networks with, e.g., 
overlaid femtocells and picocells [9], [10]. Motivated by these 
considerations, a new abstraction model is currently emerging 
and gaining popularity, according to which the positions of 
the BSs are modeled as points of a Poisson Point Process 
(PPP) and powerful tools from applied probability, such as 
stochastic geometry, are leveraged to develop tractable inte- 
grals and closed-form mathematical frameworks for important 
performance metrics (e.g., coverage and average rate) [15]- 
[17]. 

B. Stochastic Geometry based Modeling of Heterogeneous 
Cellular Networks 

The stochastic geometry based abstraction model for the 
analysis of cellular networks dates back to (at least) 1997 



IEEE TRANSACTIONS ON COMMUNICATIONS 



2 



[3], [18]. Subsequently, a similar shotgun-based, i.e., PPP- 
based, abstraction model was proposed in [19], and it was 
shown that, compared with the traditional hexagonal grid 
model, the shotgun approach provides upper performance 
bounds. More recently, the PPP model has been used for the 
analysis of spatial and opportunistic Aloha protocols [20], and 
for the characterization of the Signal-to-lnterference-Ratio 
(SIR) of (single-tier) cellular networks [21]. In spite of these 
initial and pioneering attempts of applying the PPP model 
and stochastic geometry to the analysis of cellular networks, 
only recently the random-based abstraction model for the 
positions of the BSs has received the attention it deserved. In 
particular, its emergence and widespread adoption for cellular 
networks analysis and design is mostly due to [22], where a 
comprehensive framework to compute coverage and average 
rate of single-tier deployments is provided. In [22], it is shown 
that the PPP model is as accurate as regular grid models, but it 
has the main advantage of being more analytically tractable. A 
comprehensive study based on real BSs deployments obtained 
from the open source project OpenCelllD [23] has revealed 
that the PPP model can indeed be used for accurate coverage 
analysis in major cities worldwide. Recent results about the 
validation of the PPP model for real BSs deployments are 
available in [24], where data collected from Ofcom, i.e., the 
independent regulator and competition authority in the UK, is 
used. Fueled by these encouraging results, many researchers 
are currently using the PPP-based abstraction model to study 
single- and multi-tier cellular networks, e.g., [25]-[47] and 
references therein. The PPP-based approach is also widely 
adopted for network interference modeling, e.g., [48]-[65]. 

C. Analytical Computation of the Average Rate: State-of-the- 
Art and Paper Contribution 

In this paper, we capitalize on the emerging PPP-based 
abstraction model for multi-tier cellular networks, and propose 
a new mathematical methodology to compute the downlink 
average rate over general fading channels. Recent papers 
have developed frameworks to compute the average rate for 
single-tier downlink [22], [39], multi-tier downlink [26], [43], 
downlink multi-cell coordination [36], [44], and single-tier 
uplink cellular networks [30]. All these papers use the same 
two-step methodological approach to compute the average 
rate, which was originally introduced in [20] and exploits the 
Plancherel-Parseval theorem: i) first, the Coverage Probability 
(Pcov) is computed; and ii) then, the average rate is obtained 
by integrating Pcov over the positive real axis [20, Eq. 
(2.12)]. Throughout this paper, this methodology is denoted by 
Pcov-based approach. Although this technique avoids system- 
level simulations, it requires, for general fading channels, the 
computation of a four-fold integral [22, Appendix C]. For this 
reason, many authors often limit the analysis to Rayleigh fad- 
ing channels and/or to interference-limited networks, where 
simplified frameworks can be obtained. Further details about 
the computational complexity of the Pcov-based approach are 
available in Section UlI-GI 

To overcome this limitation, we propose a new analytical 
framework which, at the same time, reduces the number of 



integrals to be computed, and, similar to the Pcov-based 
approach, is flexible enough for application to arbitrary fading 
distributions (including correlated composite channel models). 
The framework leverages the application of recent results 
on the computation of the ergodic capacity in the presence 
of interference and noise [66]. It avoids the computation of 
Pcov, and needs only the Moment Generating Function (MGF) 
of the aggregate interference at the probe mobile terminal. 
Throughout this paper, this framework is denoted by MGF- 
based approach. We show that it is applicable to multi-tier 
cellular networks with long-term averaged maximum biased- 
received-power tier association, and that either a single or a 
two-fold numerical integral need to be computed for arbitrary 
fading channels. 

D. Paper Organization 

The remainder of this paper is organized as follows. In 
Section [III the system model is described. In Section |III1 
the MGF-based approach is introduced for single-tier cellular 
networks. In Section |IV] the proposed methodology is applied 
to multi-tier cellular networks with flexible (biased) tier asso- 
ciation. In Section |V] extensive Monte Carlo simulations are 
shown to substantiate the proposed mathematical framework 
for various fading channel models and cellular networks 
deployments. Finally, Section [Vl] concludes this paper. 

II. System Model and Problem Statement 

We consider a downlink heterogeneous cellular networks 
model similar to [31], [43], and [46]. However, the following 
differences hold. In [31] and [46], the tier association policy 
is based on the instantaneous Signal-to-Interference-plus- 
Noise-Ratio (SINR). On the other hand, similar to [43], we 
consider a biased long-term averaged tier association policy, 
as described in Section Bl-BI Compared with [31], [43], and 
[46] the analytical methodology to compute the average rate is 
not based on the Pcov-based approach but on the MGF-based 
approach. 

Notation: E {•} denotes the expectation operator M.x {s) — 
E{exp(-sX)} is the MGF of random variable X. fx{-) 
denotes the Probability Density Function (PDF) of random 
variable X. T [x) = exp (— t) t^~^dt is the gamma 

function, erfc (x) = exp {— is the com- 

plementary Gauss error function. Sa,b (•) is the Lommel func- 

(&,) . 

Meijer G-function defined in [68, Sec. 2.24]. A(n,a;) — 
[x/n, [x + 1) / n, . . . , , [x + n — l)/n], with n being a positive 
integer and x a real number. The i-th entry of A (n, x) is de- 
noted by Ai {n,x). 2F1 (•, •, •, •) is the Gauss hypergeometric 
function defined in [69, Ch. 15]. iFi (•, •, •) is the confluent 
hypergeometric function defined in [69, Ch. 13]. (•!) is the fac- 
torial operator F (z, x) = /^°° exp {—t} dt is the upper- 
incomplete gamma function. 7 (z, x) — t^~^ exp {~t} dt is 
the lower-incomplete gamma function. S {■) is the Dirac delta 
function. J,^ (•) is the modified Bessel function of the first kind 
and order v defined in [69, Sec. 9.6]. ) is the Heaviside 
function, i.e., H (x) = 1 if a; > and H (x) = if a; < 0. 



tion defined in [67, Sec. 7.5.5]. 



(•) 



is the 
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pFq {■,-,■) is the generalized hypergeometric function defined 
in [70, Ch. IV]. j — \/^^ is the imaginary unit. card{-} 
denotes the cardinality of a set. 



A. Heterogeneous Cellular Networks Model 

Let us consider the PPP-based abstraction model for the 
positions of the BSs in a bi-dimensional plane. Then, a 
heterogeneous cellular deployment can be modeled as a T-tier 
network where each tier models the BSs of a particular class. 
Each class of BSs is distinguished by its spatial density (Xt 
for t= 1,2,..., T), transmit power (Pt for t = 1, 2, . . . , T), 
path-loss exponent (at > 2 for t — 1,2,..., T), biasing factor 
(Bf > for t — 1,2,...,T), and fading parameters and 
distribution. The BSs of each PPP are assumed to have the 
same transmit power, the same path-loss exponent, the same 
biasing factor, and their fading channels are independent and 
identically distributed (i.i.d.). The extension to correlated and 
identically distributed (c.i.d.) fading is discussed in Section 
IIII-FI However, for mathematical generality, we assume that 
the fading distribution of the serving (tagged) BS is different 
from the fading distribution of the intra-tier interfering BSs. 
The BSs of each tier are assumed to be spatially distributed 
according to a homogeneous PPP ($t for t ~ 1,2, ... ,T). 
The T PPPs are assumed to be spatially independent. Our 
analysis applies to a typical mobile terminal, as permissible 
in any homogeneous PPP according to the Slivnyak-Mecke's 
theorem [17, vol. 1, Theorem 1.4.5]. Without loss of generality, 
the typical Mobile Terminal (MTo) is assumed to be located 
at the origin of the bi-dimensional plane. The 6-th BS of the 
t-th tier is denoted by BSt^;,. The distance from BSt^f, to MTq 
is denoted by dt^b- The standard path-loss function I (dt^b) = 



is considered. The power channel gain of the BSt b-to- 
MTq link is denoted by gt j, = \ht^b\^, where /it f, is the related 
complex amplitude channel gain. For a fair comparison among 
fading channels with different distributions, the normalization 
constraint Ej^t.;,} = E||/if.b|^| = O = 1 is assumed for 
every b and for t — 1,2, . . . ,T. 

The frameworks developed in the present paper are ap- 
plicable to single-input-single-output transmission systems. 
In other words, BSs and MTo are equipped with a single 
transmit and receive antenna, respectively. The generalization 
of the proposed analytical methodology to more advanced 
transmission technologies is currently under investigation, but 
it is beyond the scope of the present paper The interested 
reader can, however, find preliminary results to the analysis 
of multi-antenna receivers and dual-hop relaying in [71] and 
[72], respectively. The main limitation of [71] and [72] is 
that cell association is not considered and that the distance 
from serving BS to probe mobile terminal is assumed to be 
fixed. Finally, we mention that the average rate is computed 
under the same assumptions as in [22, Sec. IV], i.e., the 
interference is treated as noise and the typical mobile terminal 
uses adaptive modulation/coding such that the Shannon bound, 
for the operating instantaneous SINR, can be achieved. 



B. Biased Long-Term Averaged Tier and BS Association 

We assume that the BSs of each tier operate in open access 
mode for MTq [8]. As a consequence, MTq is allowed to 
access to any tiers without any restrictions. In a multi-tier 
cellular networks model, both tier and BS associations have to 
be properly defined. Similar to [43, Sec. II-A], throughout this 
paper we consider a long-term averaged maximum biased- 
received-power association policy. Let dt = mm{dt.b} for 
t = 1, 2, . . . , T be the distance from MTo to the nearest BS 
of the t-th tier. Let BSt for t = 1, 2, . . . , T be the T nearest 
BSs. Then, MTq is associated (tagged) to the tier t* defined 
as follows: 

r - argmax {Ptdr"'Bt} (1) 

t=1.2,...,T 

and the tagged (serving) BS is denoted by BSf — BSq. 

In other words, MTq is connected to the BS that offers the 
highest average received power to it. Accordingly, the BSo- 
to-MTo link is the useful signal, while all the other BSs in 
every tier act as interferers. Since the positions of the BSs are 
random, the BSo-to-MTo distance is a random variable as 
well [22]. 

The biasing factor, Bt > for t = 1,2, . . . ,T, modifies the 
coverage range of each tier for a better offloading strategy. 
For example, if Bt > 1 the coverage range of the t-\h tier is 
increased. Throughout this paper, we assume, similar to [43], 
that all the BSs are fully-loaded (i.e., their queues are full 
and, thus, they have always data to transmit). The analysis 
of heterogeneous cellular networks with partially-loaded BSs 
is postponed to future research, for example either using the 
conditionally thinning approach proposed in [35] or the recent 
results in [45] and [47]. 

C. Problem Statement 

The main objective of this paper is to compute the average 
(ergodic) rate of a heterogeneous cellular network, which 
is modeled as the superposition of T independent PPPs. 
According to [22] and [43], the average rate can be written as 
follows: 



n 



T 

E 



AtTZt 



(2) 



where: i) At is the probability that MTo is associated to the 
t-th tier; and ii) TZt is the average rate of MTq conditioned on 
its association to the t-th tier For the tier association policy 
introduced in Section Hl-BI At is available in [43, Lemma 1]. 
On the other hand, TZt is defined as follows [22, Sec. IV], [43, 
Eq. (46)]: 



Ut - 



+ 00 



nt iOfd^.oiOd^ 



(3) 



where: i) dt^ is the distance of MTq from its serving BS by 
conditioning on MTq being tagged to the i-th tier; ii) fd-, g (•) 
is the PDF of the random distance dt^, which is given in [43, 
Lemma 3]: 

C T 



* I 9=1 



(4) 
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Ht = 27rAt / ^ cxp < -TT 



q=l 



(6) 



7^- / [1 - Xo (SNRy)] ^i^dy 



Gi (y) 



1 



Zi (SNRy) 2 Z/ (SNRy) 



r+oo 

/ ^^^1 exp {-ttAZ/ (SNRy) f } exp {-y^t } 
Jo 



(8) 



and iii) 7?.t (^) is the average rate of MTq conditioned on this 
terminal being tagged to the t-th tier and on c?t being equal 
to dtM = From [22, Sec. IV] and [43, Eq. (14)], TZt (•) can 
be written as follows: 



this case, ^ simplifies as follows: 



TZt (0 =IE{ln(l 



E On 1 



siNRi im 



a 

T 

9=1 f)e$5{\BSt 



o(?)} 



agg (0 



(5) 



-"5 

b 



where: i) aj^ is the noise power; ii) BSt.o {£,) is the serving 
BS at distance dt^ = £, and gt,o is the BSt_o-to-MTo power 
channel gain; and iii) /agg {0 is the aggregate interference 
conditioned on = C which is generated by all BSs except 
the serving BSt.o- From (|2]l-(|5]l, the average rate reduces 
to 7?. = Etli ^i'^h 'R-t given in ^ at the top of this page. 

The main objective of the next sections is to introduce 
a new MGF-based approach to efficiently compute TZt in 
(|6]l for arbitrary fading channels. The main contribution is 
to avoid the computational complexity of the state-of-the- 
art Pcov-based approach [20], [22], [43]. To this end, we 
introduce the simplified notation as follows, which originates 
from the assumption of identically distributed fading in each 
tier: i) ftfl (•) and Mt.o (•) are PDF and MGF of gt.o in ©, 
respectively; ii) ft.i (•) and Aitj (•) are PDF and MGF of gt^b 
in respectively; iii) Mi._^,^ (sO is the MGF of /agg (^) 
in ©; and iv) Xgj.,, (•; is the MGF of /g,agg (0 = 
Ebe<i.,{\BSto(?)} {P<i9q,bd~"''\ i.e., the per-tier aggregate 
MGF in ©. 

Finally, we mention that the average rate in (|6]l provides an 
estimate of the mean data rate over a cell that is achievable 
by a typical mobile terminal [22, Sec. IV]. This interpretation 
immediately follows from the validation procedure of the 
PPP-based abstraction model against conventional grid-based 
abstraction models, as discussed in detail in [22, Sec. V-A]. 



III. Single-Tier Cellular Networks 

To better introduce the proposed MGF-based analytical 
methodology to compute the average rate in (|6]l, we start by 
considering the single-tier reference scenario with T = 1. In 

'Throughout this paper, the serving BS is denoted by BSt,o {£.) when used 
in equations, and by BSt,o when used in the text. 



n 



2ti\ 

+ C30 







Cexp{-7rAC^}Ehn (1 



'N 



^agg 



(0 



(7) 



with /agg (0 = E6e<i.{\BSo(0} iPdbd^ Since T = 1, for 
ease of notation, in (|7]i the subscript t that denotes the tier is 
dropped. Likewise, the subscript t is dropped in /q (•), A^o (')' 
// (•), and Mi (•) as well. 

By using the MGF-based approach, an integral closed-form 
expression of (|7]i is given in Theorem [T] 

Theorem 1: Let SNR = P / ^"m '^^e Signal-to-Noise- 
Ratio (SNR), then the average rate, TZ, of a single-tier cellular 
network over generalized fading channels is given in (|8} at the 
top of this page, where: 

Zi{y)=Mi{y)+Ti{y) 

Ti (y) = r (1 - f ) E/^'mP (y) [r (2 - f + k)]"' 

Mf^ (jy)=E{5r^exp{-y5b}} 

(9) 

Proof: See Appendix H □ 

The framework in dSJ and (|9) is called MGF-based approach 
because TZ can be directly computed from the MGFs of 
useful and interference links. In fact, A^^*^^ (•) can be obtained 
from the {k + l)-th derivative of Mi {■), i.e., Mf^ (y) = 
{-d/dyf'^^ Mi{y) [73, Eq. (1.L2.9)]. In the sequel, we 
show that (•) can be explicitly computed in closed-form 

for many fading channel models. Furthermore, closed-form 
expressions of Mq (•) and Mi (•) are available in [75, Sec. 
2.2], [76, Tables II-IV], and [77, Tables II-V] for many fading 
channel models. Compared with the Pcov-based approach in, 
e.g., [20], [22], and [43], the framework in (|8]l reduces the 
number of fold integrals to be computed from four to two. 

By carefully looking at (|9]l, some important conclusions 
about the system behavior as a function of the BSs density. A, 
can be drawn, as summarized in Remark [T] 

Remark 1: Since the integrand function of Gi {■) in ((8) is 
always greater than zero, it follows that 7^ is a monotonically 
increasing function of A. Furthermore, TZ is upper-bounded as 
follows: 

1 - Mo {z) dz 



p+co 

TZ< lim 7^ (A) = 7^(^-) = / . , , , ^ , , 
~ Jo Mi (z) + Ti{z) z 

(10) 

The analytical derivation of ( fTOt is available in Appendix 
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(7rA)Zj(SNRy) \ 
2v^ 
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1/3 



27y 



Si (2/)L/2=ajv/ctD ^/(SNRy) 



'(SNRy) 



(27r)- 



a°0(7i-A)°iV2;"(SNRy) 



A(aAr,-Z/Q) 

A(aD,0) 

(11) 



nil From ([Tol l, we observe that: i) 7?.'^°°' is independent of 
the SNR = P/ajf. Thus, for very dense BSs deployments 
increasing the transmit-power does not help in increasing the 
average rate; and ii) the existence of a finite upper-bound for 
increasing A confirms that the deployment of many BSs is not 
sufficient to achieve very high data rates but more advanced 
interference management techniques seem to be needed. □ 

In the remainder of this section, we show that the two- 
fold integral in ^ can often be reduced to a single integral, 
since closed-form expressions of Gi {■) exist for many path- 
loss exponents a. Also, we show that the infinite series in 
(|9]l can be calculated for common fading distributions of the 
interference channels. 

Let us consider the computation of Gi {■) as a function of 
the path-loss exponent a. The main result is summarized in 
Corollary [T] 

Corollary 1: Let a — A, a — 6, and a/2 — a^/aD with 
aN and au being two positive integer numbers, then Gi {■) 
in ^ has closed-form expression shown in (fTTT i at the top of 
this page, where = q;/2 — L 

Proof: Equation (fTTI) follows from some notable integrals. 
More specifically: (a) from [73, Eq. (2.2. L8)]; (b) from [73, 
Eq. (2.2.L13), Eq. (2.2.L14)]; (c) from [73, Eq. (2.2.1.22)]. 
This concludes the proof. □ 

Since the case study a/2 = aN/cuD encompasses many 
scenarios of practical interest, when referring to Corollary [T] 
we will implicitly assume the closed-form expression of Gi {■) 
using the Meijer G-function. 

The single-integral expression in Theorem [T] can be effi- 
ciently computed by using the Gauss-Chebyshev quadrature 
rule, as summarized in TJemarA; |2] as follows. 

Remark 2: By using Gauss-Chebyshev integration, TZ in 
Theorem\l\can be computed as [69, Eq. (25.4.39)]: 

Ngcq 

7^~ ^ [1 _ ^0 (SNRs„)] ^7 (s„) (12) 

-I -^n 

where Wn and s„ for n = 1,2, . . . , Nqcq are weights and 
abscissas, respectively, of the quadrature rule [74, Eq. (22) 
and Eq. (23)]: 







4A^GCQ C0S2 







— tan 





^ 2n - 1 \ 


TT 


J cos 










+ 4 



(13) 



□ 



A. Computation ofTi{-) in (O/or General Fading Channels 

Theorem^and Corollary\l\need the computation of 7/ {■), 
which depends on the fading distribution of the interference 
channels. As mentioned in Theorem^ Ti (•) can, in principle, 
be computed from the derivatives of A4i (•). However, closed- 
form expressions can be obtained for many fading channel 
models by also avoiding the computation of the infinite series 
in (|9|l. Some key case studies are analyzed in Propositions\l\^ 
for Nakagami-m, Log-Normal, composite Nakagami-m and 
Log-Normal, and composite Rice and Log-Normal fading, 
respectively. 

Proposition 1: Let the interference links experience 
Nakagami-m fading. Accordingly, gt follows a Gamma 
distribution with parameters (m, fl), which we denote as 
gb ^ Gamma(m,r2) [75, Sec. 2.2.L4]. Then, 7/ (•) in ^ 
has closed-form expression as follows: 



r/(y) = m(-j 



1 - 



X 2F1 m + 1,1,2 



(14) 



Proof: See Appendix Hill □ 

Proposition 2: Let the interference links experience Log- 
Normal fading. Accordingly, gi, follows a Log-Normal distri- 
bution with parameters (in dB) (/i, cr^), which we denote as 
gb - LogN(/i,cr2) [75, Sec. 2.2.2]. Then, 7} (•) in ^ has 
closed-form expression as follows: 



Tiiy) 



r, \ -1 1 NghQ 



X exp 
X iFi ( 1,2 



|_-^Q(y2a5„+M)/l0y| 



(15) 



^ -[^q(V2(TS„+ai)/10 



where Wn and s„ for n — 1,2,..., A'^qhq are weights and 
abscissas, respectively, of the Gauss-Hermite quadrature rule 
[69, Eq. (25.4.46)]. 

Proof: See Appendix HVl □ 

Proposition 3: Let the interference links experience com- 
posite Nakagami-m and Log-Normal fading. Accordingly, gb 
follows a Gamma distribution by conditioning on its mean 
power, which, in turn, follows a Log-Normal distribution. We 
denote this distribution as gb ^ Gamma/LogN (m, /z, cr^) 
[75, Sec. 2.2.3.1]. Then, 7} (•) in ^ has closed-form expres- 
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sion as follows: 

Ti {y) « f 1 - -) 2/^ 



Nghq 
n=l 



-(m+1) 



(16) 



X 2F1 ^r7T. + l,l,2- ^,2/(2/ + m(I>„) ^ 

with (i„ = io-(v^-~"+A')/io. 

Proof : See Appendix FVl □ 
Proposition 4: Let the interference links experience com- 
posite Rice and Log-Normal fading. Accordingly, gi, fol- 
lows a non-central Chi-Square distribution by condition- 
ing on its mean power, which, in turn, follows a Log- 
Normal distribution. We denote this distribution as gi, ^ 
ChiSquare/LogN(i^,^,cr2) [73, Eq. (6)], with K being 
the Rice factor If K 7^ 0, 7/ (•) in (|9]l has closed-form 
expression given in ( fTTj i at the top of this page. If K ^ 0, 
the composite Rice and Log-Normal fading reduces to the 
composite Nakagami-m and Log-Normal fading with m = 1 
(Suzuki distribution [75, Sec. 2.2.3.2]) and Proposition^ can 
be used. 

Proof : See Appendix IVll □ 
From ( fTTT i. we observe that, unlike the other fading distribu- 
tions in Propositions\l\^ for composite Rice and Log-Normal 
fading we still need to calculate an infinite series to compute 
7/ (•) in (|9]l. The computation of the series can be avoided as 
suggested in Remark^as follows. 

Remark 3: Using the mapping between the m parameter 
of a Nakagami-OT distribution and the K factor of a Rice 
distribution [75, Eq. (2.26)], SvH can be approximated by (fTSI l 
with m = (1 + Kf ^{1 + 2K). □ 
Finally, we would like to emphasize that the fading dis- 
tributions studied in Propositions [T]-0] are just some selected 
examples, which have been chosen because they are often used 
in theoretical analysis. However, our analytical methodology 
to compute 7/ (•) in (|9|l is applicable to arbitrary fading 
distributions as described in Remark^ 

Remark 4: From (|9]l, we observe that 7/ (•) needs the 
computation of Mf\y) — E {g^'^^ exp {—ygt}} — 



exp {—yx} fg^ (x) dx. With the exception of the 



Log-Normal distribution, which is studied in Proposition |2] 
from [75, Sec. 2.2], [76, Tables II-IV], and [77, Tables II-V] 
we note that two general situations can arise: 

1) fix;y) = exp {-yx} fg^ (x) = A exp {- B (y) x} , 
where A is a constant and B {■) is a function of y. 
In other words, /(•;•) is still an exponential function 
in X. In this case, closed-form expressions of 7} (•) 



can be obtained by using the same development as in 
Proposition [T] for Nakagami-m fading. 



2) fsA^) 



LjX 



Dx 



, where C, D, and 



V are constants. In other words, the distribution of the 
power channel gain gi, can be cast in terms of a Meijer 
G-function. Accordingly, Ad^'^'' (•) can be computed in 
closed-form as another Meijer G-function by using the 
Mellin-Bames theorem and the notable integral in [68, 
Eq. (2.24.3.1)]. In general, in this case it is not possible 
to avoid the computation of the infinite series in (|9]l. □ 

B. Efficient Computation of the Meijer G-Function in ( 1771 ) 

The computation of the average rate in (|7]l by using Theorem 
[U and Corollary [T| needs, in general, the calculation of the 
Meijer G-function in (fTTT i. This special function is commonly 
used in wireless communications theory, e.g., [76], [77], [79]- 
[81], and it is available in several standard mathematical 
software packages. Thus, in general, its computation can be 
performed very efficiently. However, in (|8) the Meijer G- 
function must be calculated for all positive real values, and 
it is known that the numerical complexity and the numerical 
accuracy of common algorithms to compute the Meijer G- 
function increases and decreases, respectively, for small values 
of its argument, i.e., for y — in ( fTTT i. see, e.g., [82] and [83]. 
In order to provide a framework that is general and accurate 
but also simple and stable to compute. Corollary |2] provides 
a numerically efficient and stable solution to compute Qi {■) 
in ( fTTT i. which exploits an asymptotic expansion of the Meijer 
G-function for large values of its argument. 

Corollary 2: Let a/2 — aM/ao with and being 
two positive integer numbers, then Gi { ) in (O can be effi- 
ciently computed as shown in ( fTsT i at the top of the next page, 
where e is a small positive constant. 

Proof: See Appendix [Vll] □ 
The rationale behind Corollary^is to avoid the calculation 
of the Meijer G-function for small values of its argument, and 
to replace the Meijer G-function with an accurate, simple to 
compute, and numerically stable expansion formula. In other 
words, Qj (•) is computed by using U (•), i.e., the exact formula 
in Corollary |2] as long as the argument of the Meijer G- 
function is no smaller than e. On the other hand, when this 
occurs the asymptotical expansion ^(asymptote) ^.-j used, 
which is simple and fast to be computed. This "adaptive" 
approach allows us to keep the desired accuracy without 
increasing the numerical complexity and without incurring in 
numerical instabilities. The key parameter for the efficient 
computation of Gi {■) in ( fTsT i is e, which depends on the 
mathematical software package being used to compute the 
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Meijer G-function. In practice, e is the smallest value of 
the argument of the Meijer G-function for which it can 
be efficiently computed. If e = 0, Corollary |2] reduces to 
Corollary [T] 

Finally, we close this section with Remark |5] and Remark 

m 

Remark 5: From ( fTSI l, we note that a very computation- 
ally efficient framework, which is accurate for sufficiently 
small values of e, can be obtained by simply neglecting 
^(asymptote) (^.-j^ shown in ([Hll at the top of this page □. 

Remark 6: The integral in Gi {■) belongs to the so-called 
"Weibull-type" integrals, since it coincides with the integral 
to be computed to obtain the MGF of the Weibull distribution 
[84, Eq. (2)] and [85, Eq. (2)]. Similar to Corollary \l\ it can 
be computed in terms of Meijer G-function [84] and [86] 
or in terms of generalized hypergeometric function [85]. Fur- 
thermore, various closed-form approximations are available 
in the literature, such as [87]-[90] and references therein. The 
interested reader can consult these papers and the references 
therein to identify alternative ways of computing Gi {■) that 
avoid special functions, such as the Meijer G-function. On the 
other hand, to the best of the authors knowledge, the approach 
proposed in Corollary^is not available in the literature. □. 

C. Interference-Limited Scenario 

In many practical situations of interest, the background 
noise is often negUgible compared to the aggregate interfer- 
ence [22] and [43]. In this case, r/!eorem[T]simplifies as shown 
in Corollary^ 

Corollary 3: Let cr^ — 0, then the average rate, TZ, in ^ 
simplifies as follows: 

Jo Mi{z)+Ti{z) z ^^^^ 

Proof: By using the change of variable SNRy = z in ([8]l, 
it follows that Gi (SNR'^z) = 1/Z/ (z) since SNR ^ oo if 
(jj^ = 0. This concludes the proof. □ 



From (|20l i, interesting considerations about the average rate 
can be made, as summarized in Remark [T] 

Remark 7: Similar to [22, Eq. (8)], (|20] | confirms that 
for interference-limited cellular networks the average rate is 
independent of the density of BSs as well as of the transmit- 
power Thus, increasing either the BSs density or the transmit- 
power are not effective solutions to increase the average 
rate. More advanced interference management mechanisms 
are needed. Furthermore, our framework shows that these two 
trends hold regardless of the fading channel model, and that 
they seem to be mainly related to the PPP spatial model of 
the BSs. Finally, by comparing ( |20l i with ( fTOl i we observe that 
TZ < 7?.(-^^) = 7^(SNRoo) xhis implies that the average rate 
of a cellular system with unbounded BSs density and finite 
transmit-power is the same as the average rate of a cellular 
system with unbounded transmit-power and finite BSs density. 
□ 

D. High-SNR Scenario 

In Corollary [3] we have studied the average rate in the 
absence of background noise. In Corollary ID we study the 
scenario with small but non-zero noise, i.e., the high-SNR 
setup. 

Corollary 4: As a function of the SNR, the average rate, 
TZ, in (|7]l is upper- and lower-bounded as follows: 

^(SNR»l) ^ 7^(SNR^) _ (^;,)-"/2 y(i + ^)J— 

^ ' V 2/ SNR 



Jo [Mj{z)+Ti{z)]'+^^/'^ 

< 7^(SNR) < 7^(^^^~) 

Proof : Equation dSTT i immediately follows, with the same 
analytical steps, from (|39] l in Appendix |ll] and by using the 
identity {a/2) T (a/2) = T (1 + a/2) in gOll. This concludes 
the proof. □ 

From (l2Tl i. interesting considerations about the average rate 
can be made, as summarized in Remark [8j 
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Remark 8: By direct inspection of ( 1211 1. we observe that: i) 
the lower-bound, 7^(snr»i)^ jj^g high-SNR approximation 
of the average rate since limsNR-^+oo 7^(SNR»l) ^ 7^(snr^). 
ii) the average rate increases with the SNR by approaching the 
upper-bound 7^(snRoo) y^jfjj linear convergence rate; iii) the 
average rate increases with the BSs density by approaching 
the upper-bound TZ(^^^^) — 7^(^=0) with (a/2)-order con- 
vergence rate; and iv) the larger the path-loss exponent, a, 
the faster the convergence speed to TZ(^^^<=^'> — 'R,'y^°^) as a 
function of A. □ 

E. Frequency Reuse 

In this section, we study the impact of frequency reuse on 
the average rate. In particular, we consider a cellular network 
with Fb > 1 frequency bands. The setup with Fb = 1 
corresponds to the universal frequency reuse case studied 
in Theorem [T] Also, similar to [22], we assume that each 
interfering BS picks at random one of the Fb frequency bands 
when transmitting. The average rate is given in Corollary^ 

Corollary 5: The average rate, TZ'^-^^K of a single-tier 
cellular network with Fb > 1 available frequency bands 
and random frequency reuse coincides with TZ in (|8]l, ( fTOl l, 
(EOll, and (EB by replacing Zj {■) in Q with zj^^' (z) ^ 
{Fb - 1) + Mi {z) + Ti (z), and A (when available) with 

XiFB)^X/FB. 

Proof : The proof follows by taking into account that for 
Fb > 1: i) the average rate in (|5) becomes TZ^^^^ (^) = 
(l/i^s) E {In (1 + SINK (C))}; ii) on average, the interference 
originates from a PPP with BSs density equal to A'^^' — 
\/Fb', and iii) the tier and BS association PDF is independent 
of Fb [22, Sec. VI-A]. Accordingly, the proof proceeds along 
the same lines as Theorem^ This concludes the proof. □ 

From Corollary |5] interesting considerations about the av- 
erage rate can be made, as given in Remark |9] 

Remark 9: By direct inspection of TZ^^"^ in Corollary^ it 
follows that the average rate is maximized for Fb — 1, i-e., 
for universal frequency reuse. Also, for densely deployed BSs 
([Tol l and for interference-limited cellular networks (E&i . the 
average rate linearly decreases with the number of available 
frequency bands Fb- These trends are in agreement with [22, 
Sec. VI-B] and hold for general fading channel models. □ 

F. Correlated Log-Normal Shadowing 

The MGF-based approach introduced so far can be applied 
to a wide variety of channel conditions, notably composite 
fading channels that account for Log-Normal shadowing (see 
Proposition^wd Proposition^. The average rate in Theorem 
[T] is applicable, however, only to i.i.d. fast-fading and Log- 
Normal shadowing. It is well-known, on the other hand, 
that shadowing correlation severely affects the performance 
of cellular networks [93]-[97]. In this section, we provide a 
simple methodology to extend the framework in Theorem [T] 
to c.i.d. fading channels. The reason of restricting the anal- 
ysis to equi-correlated fading originates from the stochastic 
geometry approach for other-cell interference modeling used 
in the present paper, which is applicable only to identically 



distributed fading. The methodology used to obtain equi- 
correlated Log-Normal random variables exploits the Owen 
and Steck method for the generation of equi-correlated mul- 
tivariate Normal distributions [98]. 

As an illustrative example, let the generic downlink chan- 
nel experience composite Nakagami-w fast-fading and Log- 
Normal shadowing, as described in Proposition |3] The pro- 
posed methodology is readily applicable to other fading 
channel models with correlated Log-Normal shadowing, as 
well as to multi-tier cellular networks by applying the same 
methodology to the framework discussed in Section |IV] More 
specifically, we assume i.i.d. fast-fading and c.i.d. shadowing. 
Accordingly, card{'l>} channel power gains gi, for 6 e $ 
with correlation coefficient p and parameters (m, /i,CT^) can 
be obtained as follows [98]: 

Step 1 : Generate card {$} equi-correlated Normal random 
vai-iables as — a^fpS-l-ay/l — pSb+fi for 6 S $, 
where S and Sb for 6 S $ are a set of i.i.d. Normal 
random variables with zero mean and unit variance. 
The card {$} random variables Xb have mean p and 
variance a^, for b E ^, regardless of the correlation 
coefficient p. 

Step 2 : Convert the set of card {$} equi-correlated Nor- 
mal random variables into a set of card{$} equi- 
correlated Log-Normal random variables as Yb = 

Step 3 : Generate card{<I>} independent Gamma random 
variables gb with fading severity m and mean value 
Yb for 6 e $. 

From the above generation mechanism of c.i.d. com- 
posite Nakagami-OT and Log-Normal fading, it follows 
that the card{$} random variables gb are, by condition- 
ing upon the random variable S, i.i.d. with parameters 
(to, p + (T^S, (1 — p)). Accordingly, we propose the fol- 
lowing approach to compute the downlink average rate of 
cellular networks: 

Step 1 : The framework in Theorem [T] is applied by con- 
ditioning upon the random variable S and by sub- 
stituting p i—i' p + (Ty/pS and (T I—)- cr-y'l — p. The 
resulting average rate is denoted by TZ (5). 

Step 2 : The conditioning upon the standard Normal random 
variable S is removed by averaging over its PDF 
/s(a;) = (l/V2^)exp{-a;72}. 

In formulas, the downlink average rate over c.i.d. fading 
channels can be computed as follows: 

r+°° (a) 1 '^^"^ 

TZ = J TZ{x)fs (x) dxK --j= ^^),,7^ \^s^ (22) 

where (a) is obtained by applying Gauss-Hermite quadratures, 
and TZ {x) is the average rate in Theorem [T] with fading 
parameters to (x) — m, p{x) = p + <jy/px, and a (x) — 

In summary, the rationale behind the proposed approach to 
deal with shadowing correlation consists in: i) first, generating 
a set of correlated Log-Normal random variables that are 
conditionally independent and, thus, applying the framework 
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for independent shadowing; and ii) then, removing the con- 
ditioning via a single numerical integration. Accordingly, 
shadowing correlation can be taken into account with only 
a single extra numerical integral, which can be efficiently 
computed using Gauss-Hermite quadratures as shown in (l22t . 

G. Pcov- vs. MGF-based Approach: A Comparison 

In Section II-CI we have stated that both Pcov- and MGF- 
based approaches can be applied to general fading distribu- 
tions. However, Pcov- and MGF-based approaches need, in 
general, the computation of a four- and a two-fold numerical 
integral, respectively. In both cases, the integrals may involve 
the computation of special functions, which, however, are 
efficiently implemented in commercially available software 
packages. Due to the reduction of the number of fold integrals 
to be computed, the MGF-based approach is expected to be 
more computationally efficient. The aim of this section is 
to better compare strengths, weaknesses, and computational 
complexity of these two approaches. 

To better conduct this comparison, we summarize in (|23), 
shown at the top of this page, the four-fold integral expression 
of the average rate that is obtained from the Pcov-based 
approach. More specifically, (|23) is obtained from [22] as 
follows: (a) originates from [22, Appendix C]; (b) originates 
from [22, Appendix B]; and (c) originates from [22, Theorem 
4, Eq. (4)]. For consistency and ease of comparison, the same 
notation as for the MGF-based approach is used. 

By comparing the MGF-based approach in dH) with the 
Pcov-based approach in ( l23T l the following comments can be 
made: 

• Both approaches may need the computation of some 
special functions. More specifically, the MGF-based ap- 
proach involves the computation of hypergeometric func- 
tions in Ti {■), and the Pcov-based approach involves 
the computation of the incomplete Gamma function in 

• Both approaches may need to use Gauss-Hermite quadra- 
tures to compute /o (•), A^o (■)> // (')' ™d Mi (•) for 
composite channel models. This need originates from the 
analytical intractability of Log-Normal shadowing and it 
is independent of either the Pcov- or the MGF-based 
approach being used. 

• By using Corollary [T] the two-fold integral in (|8]l may 
be reduced to a single-integral for some path-loss ex- 
ponents. Likewise, by using the Mellin-Barnes theorem 
in [68, Eq. (2.24.2.1)] and the Meijer G-function rep- 
resentation of the upper-incomplete Gamma function in 
[68, Ea. (8.4.16.2)], a closed-form expression of (•, •) 



in ( l23T l may be obtained. As a consequence, the MGF- 
based approach reduces the number of fold integrals to 
be computed and avoids the computation of complex 
integrals. 

• In interference-limited scenarios, the MGF-based ap- 
proach in Corollary [3] offers a significant reduction of 
the computational complexity and the average rate can 
be calculated from the simple single integral in (|20] |. 
On the other hand, the computational complexity of 
the Pcov-based approach is not significantly affected in 
this scenario. In fact, the only simplification in ( |23] l is 
exp {— 27rcr^js} = 1 in Ti (•, •), which does not lead to 
further reduction of the number of fold integrals to be 
computed. 

• The desired form of the average rate offered by the MGF- 
based approach leads to simple and intuitive understand- 
ing of the performance of cellular networks for a variety 
of special operating scenarios, such as dense cellular net- 
works deployments {Remark [T]), interference-dominated 
environments {Remark^wA Remark^, frequency reuse 
strategies {Remark |9]l. On the other hand, little insight 
can be gained from ( |23] l for general fading distributions. 
However, (|23] | can be significantly simplified for Rayleigh 
fading channels and interesting design guidehnes can be 
inferred from it [22]. 

The considerations above originate from the direct inspec- 
tion of (O and ( |23] ), and provide a quahtative comparison 
of the reduction of computational complexity that can be 
expected by using the MGF- instead of the Pcov-based 
approach. To better understand the advantages of the MGF- 
based approach, we have also conducted some numerical tests 
with the goal of providing a more quantitative assessment of 
the computational complexity. The conventional approach that 
is often used to conduct these tests it to consider a case study 
for which, with further analytical manipulations, the integral 
expressions in ([8]l and ( |23] l can be simplified or even computed 
in closed-form, and to compare their accuracy and computa- 
tional time without applying any mathematical manipulations. 
Following this line of thought, we have considered Rayleigh 
fading as a benchmark and have implemented in Mathematica 
the formulas in (|8]l and ( |23] l as they appear in the present paper 
In fact, simple closed-form expressions for Rayleigh fading 
are available in [22]. As far as the MGF-based approach is 
concerned, the outer integral in (|8) is computed using (ITZt 
with A^GCQ = 2000. The high value of iVccQ is chosen as a 
worst case setup for the MGF-based approach. Various com- 
binations of path-loss exponents, a = {2.05, 2.2, 2.5, 3, 4, 5}, 
and densities of BSs, A = {lO~^ lO""*, 10"^ 10-^}, have 
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been considered. The chosen path-loss exponents cover typical 
propagation environments for cellular applications [99, Table 
2.2], [100, Ch. 2, Sec. 5], and the chosen densities of BSs 
cover sparse, normal, and dense cellular deployments [22], 
[42]. The tests have been executed in a laptop computer In 
all the analyzed scenarios, the MGF-based approach in (|8) 
has been able to provide accurate estimates of the average 
rate in less than five/six seconds for each SNR point to be 
computed. On the other hand, the Pcov-based approach in 
( |23] | has not been able to provide any numerical estimates after 
five minutes of computation. In interference-limited scenarios, 
i.e., (7^ = 0, the computational complexity of the MGF- 
based approach is further reduced, while the computational 
complexity of the Pcov-based approach is not affected. These 
outcomes confirm the advantages of the proposed MGF-based 
approach for analysis and design of cellular networks. 

IV. Multi-Tier Cellular Networks 

In this section, we extend the analytical framework to 
generic multi-tier cellular networks. The analytical develop- 
ment is, in many ways, similar to Section HUl Thus, only the 
most important analytical details are reported in what follows. 
The departing point is (|6]l and the main result is summarized 
in Theorem^ 

Theorem 2: Let for t = 1, 2, . . . , T in ©. Let SNR* = 
Pt/<j'jq be the SNR of the i-th tier, an explicit closed-form 
expression of TZt for arbitrary fading channels is given in 
at the top of this page, where: 




SNR,y 



SNR,y 



%j (y) - r 



1-A 



+ CX3 

X Vy 



fc=0 



'^'M^^} (y) 



r 2- 



(25) 



(26) 



and M^} is) = E j^+i exp {-sg,,b}}. 

Proof : The proof follows by using the same steps as in 
Appendix U and by taking into account that, thanks to the 



spatial and channel independence of the PPPs, the MGF 
of the aggregate interference for the generic <-th tier is 

-^/a,, = nr=i>'9,/a,, (s;^?) where M^j^^^ is 
given in ( l36l ) and can be computed in closed-form from Xq, 

Mqj (•), %,i (■), and = iPq/Ptf"^ {BjBtf-^- e*/"' 
[43, Eq. (42)]. This concludes the proof. □ 

Theorem |2] provides a very general expression for the 
average rate of multi-tier cellular networks that, in general, 
needs the computation of a two-fold numerical integral but 
is applicable to tiers having different path-loss exponents and 
fading distributions. Furthermore, closed-form expressions for 
Tqj (•) in ( l26b can be obtained from Propositions^^^ similar 
to the single-tier case. The extension to correlated Log- 
Normal shadowing follows immediately from Section IIII-FI 

Even though general, r/ieoremll] provides a framework that 
is less analytically tractable than Theorem [T] and Corollary 
[U A simpler and more insightful analytical framework can be 
obtained by assuming that all the tiers have the same path-loss 
exponent, i.e., at = a for t = 1,2, ... ,T while still keeping 
the assumption that the per-tier fading distribution is different. 
The related framework is given in Corollary |6l 

Corollary 6: Let at = a for t = 1,2, ... ,T, then TZt in 
(|6]l can be explicitly computed as shown in dZTI l at the top of 
this page, where: 




] SNR,yj 

SNR^y^J 

(28) 

Proof: The proof follows directly from (l24t and from some 
algebraic manipulations similar to Appendix |I] This concludes 
the proof. □ 

From Corollary |6] the following important remark can be 
made. 

Remark 10: By comparing Corollary^and Corollary\l\ we 
observe that the two formulas have the same structure. More 
specifically, TZt in ( l27b can be obtained from 72, in (O by 
simply replacing Zj (•) with zj*"* (•). Furthermore, Theorem^ 
reduces to Corollary^ as expected, for T = 1. □ 
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Remark \TU\ allows us to easily generalize many important 
results obtained for the single-tier setup to the multi-tier case. 
In particular. Corollary |7] generalizes Corollary |5] for multi- 
tier cellular networks with random frequency reuse; Corollary 
[8] generalizes Remark [T] by investigating the impact of dense 
BSs deployments; and Corollary |9] generalizes Corollary [3] 
and Corollary |4] by studying interference-limited multi-tier 
cellular systems and high-SNR operating conditions. 

Corollary 7: The average rate of the t-th tier of a multi- 
tier cellular network with at = a and with F^*-* > 1 available 
frequency bands and random frequency reuse for every tier t = 
1, 2, . . . , T can be obtained from ( |27| i by replacing zj*''''' (•) 
in (|28]l with: 

il''" (y) = {"L" - 1) 

(29) 

for t,q = l,2,...,r, and Af with a} ^ = Xt/Pj^^ for 
1,2,...,T. 

Proof: It follows from Corollary |5] and Remark [TO] This 
concludes the proof. □ 

Corollary 8: Let us consider a multi-tier cellular network 
with at = a and with > 1 available frequency bands 
and random frequency reuse for every tier t = 1,2, ... ,T. 
Also, let Ai = KtX and SNR* = XtSNR for t = 1, 2, . . . , T. 
If A — > +00, then the average rate in Corollary [T] is upper- 
bounded as shown in (l30l l at the top of this page. 

Proof : The proof follows immediately from Remark [T] and 
Remark\lU\ This concludes the proof. □ 

Corollary 9: Let us consider a multi-tier cellular network 
with at — a and with Fg^ > 1 available frequency bands 
and random frequency reuse for every tier t ~ 1,2, ... ,T. 
Also, let At = KtX and SNRt = XtSNR for t = 1, 2, . . . , T. 
If ajf = (i.e., SNR — > +00), which implies that the system 
is interference-Umited, then the average rate in Coro/Zary [T] is 



upper-bounded as iZt < = iz['^°°\ where 

is given in dSOl l. Furthermore, if aj^ is small but non-zero, 
i.e., SNR ^ 1, then the average rate in Corollary ^can be 
approximated as shown in dSTT l at the top of this page. 

Proof: The proof follows from Corollary [3] Corollary |4] 
and Remark\lU\ This concludes the proof. □ 

V. Numerical and Simulation Results 

In this section, we show some numerical examples in order 
to verify the accuracy of the proposed analytical methodol- 
ogy against Monte Carlo simulations, as well as to show 
the impact of different fading parameters and distributions 
on the average rate. For a fair comparison among different 
fading distributions, the mean square value of each fading 
distribution is normalized and set equal to one. This implies 
n = 1 for Rayleigh and Nakagami-m distributions, and 
/i — — In (10) cr^/20 for Log-Normal, composite Nakagami- 
m and Log-Normal, and composite Rice and Log-Normal 
distributions. Furthermore, useful and interference links are 
assumed to have the same fading distribution. The analytical 
framework is implemented as described in the captions of 
each figure. As far as the composite Rice and Log-Normal 
fading model is concerned, both frameworks in Proposition |4] 
and Remark |3] are implemented. We have verified that both 
frameworks provide the same accuracy. Thus, the application 
of Remark^is recommended since it is simpler to compute. 

a) Monte Carlo Simulations: As far as Monte Carlo sim- 
ulations are concerned, we have used the following method- 
ology [46, Appendix F]. 

Step 1 : A finite circular area of (normalized) radius Ra 
around the origin, i.e., where the probe mobile ter- 
minal is located, is considered. The radius is chosen 
sufficiently large to minimize the error commit- 
ted in simulating the infinite bi-dimensional plane. 
In the considered setup, the radius Ra is cho- 
sen such that Amin^^A — 100, where Amin = 
min {Ai, A2, . . . , At}- For example, Ra = 100 if 
Amin = 10"^ and A^in = 10"^ Ra = 1000 if 
Amin 10-"*, and Ra = 10000 if A^in = 10'^- 
Step 2 : For each tier, the number of BSs is generated 
following a Poisson distribution with density At and 
area tiR\. 
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Step 3 : The BSs of each tier are distributed following a 
uniform distribution over the circular region of area 
nR\. 

Step 4 : Independent channel gains are generated for each 

BS of every tier 
Step 5 : The tier and BS association policy described in 
Section III-BI is applied, and useful and interference 
Unks are identified. 
Step 6 : Given the associated tier and its tagged BS, the 

SINR is computed as shown in (|5]l. 
Step 7 : The rate of the generic Monte Carlo trial is com- 
puted as 7^mc = (l/i^4**^)ln(l + SINRt.), where 
t* is the tagged tier 
Step 8 : Finally, the average rate is computed by repeating 
Step 1-Step 7 for A^nic times and eventually calculat- 
ing TZ = (l/iVmc)I]mc=i'^nir. In our simulations. 



a = 5 



we have considered = 10^. 



In Section IIII-GI we have compared the computational 
complexity of Pcov- and MGF-based approaches, and we 
have shown that the proposed analytical methodology turns 
out to be more computational efficient for general fading dis- 
tributions. As far as the computational complexity comparison 
with Monte Carlo simulations is concerned, our experiments 
have revealed that each simulation curve shown in this section 
can be obtained in a computation time of the order of a few 
minutes (five to ten minutes depending on the setup) by using 
the MGF-based approach. On the other hand, the same curve 
can be obtained in tens of hours (ten to sixty depending on the 
setup) of computation time by using Monte Carlo simulations. 
In addition to the longer simulation time and to the more 
resources for the computation, it is important to mention that 
Monte Carlo simulations tend to be less accurate for: i) sparse 
cellular networks; ii) low path-loss exponents; and iii) high- 
SNR. The reason is that in these operating scenarios Ra and 
A^mc must be increased in order to account for the interfering 
BSs that are far from the probe mobile terminal, and which, 
in these cases, can no longer be neglected. 

b) Framework Validation for Single-Tier Cellular Net- 
works: In Figs. [T]-[5] the average rate of Rayleigh, Nakagami- 
m, Log-Normal, composite Nakagami-m and Log-Normal, 
and composite Rice and Log-Normal fading is shown, respec- 
tively, for a single-tier cellular network. Overall, we observe 
a very good accuracy of the proposed MGF-based approach. 
Furthermore, we observe, as expected, that the average rate: 
i) increases with the BSs density; ii) depends on the fading 
distribution; and iii) increases with the path-loss exponent. 

In Figs. [T}{5] we have considered dense cellular networks 
(A > 0.1) and large path-loss exponents (a > 4). The 
reason of this choice is mainly due to the long time needed 
to obtained Monte Carlo simulations for less dense cellular 
networks and for smaller path-loss exponents. However, it is 
important to verify the accuracy of the proposed MGF-based 
approach for more practical densities of BSs and for a wider 
range of path-loss exponents. In general, practical densities 
for macro BSs deployments are of the order of A w 10^^ 
[22], [42], [43]. Thus, to test numerical accuracy and stability 
of the MGF-based approach, we consider densities of BSs 
in the set A = {lO"^ 10"'', 10"^ lO^^}, in order to study 
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Fig. 1. Average rate of a single-tier cellular network over Rayleigh fading 
(Fg = 1). Markers show Monte Carlo simulations. Solid lines show the 
analytical framework, which is computed by using Theorem [T] Corollary [2] 
with e = 0.05, ProposiYi'on [T] with (m = 1, H = 1), and Remark [2] with 
^GCQ = 2000. Furthermore, Mo (■) = Mi (■) are obtained from [75, Eq. 
(2.8)]. 
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Fig. 2. Average rate of a single-tier cellular network over Nakagami-m 
fading (Fg = 1). Markers show Monte Carlo simulations. Solid lines show 
the analytical framework, which is computed by using Theorem^ Corollary 
|2]with £ = 0.05, Proposition^m\h (m = 2.5, Q. = 1), and Remcirt [2] with 
^GCQ = 2000. Furthermore, Mo (■) = Mi (■) are obtained from [75, Eq. 
(2.22)]. 



sparse, medium, and dense deployments. As far as the path- 
loss exponent is concerned, we consider values in the set a — 
{2.05,2.2,2.4,3,4,5} > 2, which cover typical propagation 
environments for cellular applications [99, Table 2.2], [100, 
Ch. 2, Sec. 5]. The results of this study are shown in Figs. 

In Figs. |6] and |2l we compare Monte Carlo simulations 
with the MGF-based approach. In particular, numerical results 
are obtained by using both the (exact) two-fold integral in 
77ieorem[T]and the (approximated) single-integral in Corollary 
|2] in order to test complexity and accuracy of both meth- 
ods. As far as Monte Carlo simulations are concerned, it 
is worth mentioning that only some SNR points (markers) 
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Fig. 3. Average rate of a single-tier cellular network over Log-Normal 
fading (Fg = 1). Markers show Monte Carlo simulations. Solid lines show 
the analytical framework, which is computed by using Theorem^ Corollary^ 
with £ = 0.05, Pmposition^witii (cr = 6 dB, = - In (10) 0-2/20 dB) and 
Nqhq = 5, and Remar*: |2] with Nqcq = 2000. Furthermore, Mo (■) = 
Mi (■) are obtained from [75, Eq. (2.54)] with A^ghq = 5. 
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Fig. 4. Average rate of a single-tier cellular network over composite 
Nakagami-m and Log-Normal fading (Fg = 1). Markers show Monte Carlo 
simulations. Solid lines show the analytical framework, which is computed by 
using Theorem^ CowUary^with e = 0.05, Proposition^with (m = 2.5, 
o- = 6dB, ^ = - In (10) 0-2/20 dB) and Nq^iq = 5, and Remart |2] with 
^GCQ = 2000. Furthermore, Mo { ) = Mi {■) are obtained from [75, Eq. 
(2.58)] with ATghq = 5. 



are shown in the figures. The missing SNR points are not 
shown because of the long simulation time and the need to 
consider very large simulation areas to get accurate estimates 
of the average rate. The simulation time increases, in general, 
for more sparse cellular networks and for smaller path-loss 
exponents. The SNR points shown in the figures are those for 
which accurate estimates can be obtained with the simulation 
setup described above. The numerical examples confirm the 
very good accuracy and the numerical stability of the MGF- 
based approach for all cellular network setups. By comparing 
the curves obtained using Theorem [T] and Corollary |2] we 
notice that the latter is very accurate except in the transition 



Fig. 5. Average rate of a single-tier cellular network over composite Rice and 
Log-Normal fading (Fg = 1). Markers show Monte Carlo simulations. Solid 
lines show the analytical framework, which is computed by using Theorem 
[T] Coro//an'|2]with e = 0.05, ProposiYion g] with (K = 10, o = 6 dB, = 
-ln(10)o-2/20dB) and ArQHQ = 5, and Remar/t |2] with Nqcq = 2000. 
In Proposition^ the series in (TT) is truncated to the first 100 terms. The 
application of Remark: [5] provides the same result and accuracy, but with less 
computational complexity. Furthermore, Mo {■) = Mi (■) are obtained from 
l|46) with A^GHQ = 5. 



a = 4 a = 5 
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Fig. 6. Average rate of a single-tier cellular network over composite 
Nakagami-m and Log-Normal fading (Fg = 1). Markers show Monte Carlo 
simulations. Solid lines show the analytical framework, which is computed by 
using Theorem^ Coro//ary |2] with e = 0.05, PropoiiYion [5] with (m = 2.5, 
o- = 6dB, /i = -In (10) 0^/20 dB) and A^ghq = 5, and Remar/t |2] with 
^GCQ = 2000. Furthermore, Mo { ) = Mi (•) are obtained from [75, Eq. 
(2.58)] with A^GHQ = 5- The black dashed lines are obtained by directly 
computing the two-fold integral in r/ieorem [T] without using the Meijer G- 
function in Corollary^ 

(corner) region from noise- to interference-limited operating 
conditions, where the average rate reaches the asymptote 
calculated in Corollary |3] (further comments are available 
below where the high-SNR scenario is discussed). The reason 
of this (in practice negligible) numerical inaccuracy originates 
from the non-smooth transition introduced by the adoption of 
the Heaviside function in Corollary |3] 

In Fig. [H we compare the average rate as a function of the 
path-loss exponent. Only numerical results obtained from the 
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Fig. 7. Average rate of a single-tier cellular network over composite 
Nakagami-m and Log-Normal fading (Fg = 1). Markers show Monte Carlo 
simulations. Solid lines show the analytical framework, which is computed by 
using Theorem\T\ Coro/Zor)' [2] with e = 0.05, Pwposition^with (m = 2.5, 
cr = 6 dB, At = - In (10) 0-2/20 dB) and AfQHQ = 5, and Remar/t |2] with 
^GCQ = 2000. Furthermore, A4o {■) = M.i {•) are obtained from [75, Eq. 
(2.58)] with Nqy{q = 5. The black dashed lines are obtained by directly 
computing the two-fold integral in T/ieorem [T] without using the Meijer G- 
function in Corollary^ 
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Fig. 8. Average rate of a single-tier cellular network over composite 
Nakagami-m and Log-Normal fading (Fb = 1). The curves are obtained 
by using the analytical framework and Monte Carlo simulations are not 
shown. More specifically, the curves are obtained by computing the two- 
fold integral in Theorem^ PropaviYion [3] with (m = 2.5, a = 6dB, /i = 
- In (10) (T2/20dB) and AfcHQ = 5, and Remar/t |2] with A^gcq = 2000. 
Furthermore, A4o { ) = A4j (■) are obtained from [75, Eq. (2.58)] with 
A^GHQ = 5. 



MGF-based approach are shown in this figure. The curves are 
obtained by using the two-fold integral in Theorem [T] As far 
as the application of Corollary |2] is concerned, its accuracy 
for a > 2.4 is shown in Figs. |6]and|7] For a < 2.4, it is less 
practical to use the Meijer G-function in Corollary |2] since 
for such values of a we would have a^v ^ 1 and ajj ^ 1, 
and, thus, computation time and numerical accuracy would 
highly depend on the specific implementation of the IVIeijer G- 
function. On the other hand, r/!eorem[T]provides very accurate 



Fig. 9. Average rate of a single-tier cellular network over Nakagami-m 
fading {Fb = 1, SNR = 10 dB). Markers show Monte Carlo simulations. 
Solid lines show the analytical framework, which is computed by using 
Theorem [T] Corollary |2] with e = 0.05, Proposition [T] with £7 = 1, and 
Remark^viilh Nqcq = 2000. Furthermore, Mo ( ) = Mj (■) are obtained 
from [75, Eq. (2.22)]. 



estimates in a few seconds (for each SNR point), as discussed 
in Section IIII-GI Figure [8] clearly shows that Theorem [T] 
provides reliable numerical estimates for the considered set 
of path-loss exponents, including a « 2. As far as the 
performance trend is concerned. Fig. |8] shows a very different 
behavior for dense (A = 10~^) and medium/sparse (A = 10^"*) 
cellular networks. In dense cellular networks, the higher the 
path-loss exponent the better the average rate regardless of the 
operating SNR. On the other hand, in medium/sparse cellular 
networks two SNR regions can be identified: i) for low-SNR 
(noise-limited regime), the lower the path-loss exponent the 
higher the average rate. This is due to that fact that the useful 
signal undergoes a lower attenuation and that the aggregate 
interference is negligible compared to the additive noise; and 
ii) for high-SNR (interference-limited regime), the higher the 
path-loss exponent the higher the average rate. This is due to 
the fact that the additive noise is negligible compared to the 
aggregate interference and that the interfering BSs undergo 
a larger attenuation, which has a more pronounced effect on 
the average rate than the larger attenuation undergone by the 
useful signal. These results are in agreement with intuition and 
confirm the usefulness of the proposed MGF-based approach 
for cellular networks analysis and design. 

For ease of comparison with Monte Carlo simulations, in the 
following only large path-loss exponents and dense cellular 
networks deployments are considered. 

c) Impact of Fading Model and Fading Parameters: In 
Figs. I9UT2I the average rate of Nakagami-m, Log-Normal, 
composite Nakagami-m and Log-Normal, and composite Rice 
and Log-Normal fading is shown, respectively, for a single- 
tier cellular network and for different choices of the fading 
parameters. Also in this case, the framework provides very 
accurate estimates. More specifically, we observe that: i) the 
average rate is slightly sensitive to m and K fading parameters 
of Nakagami-m, composite Nakagami-m and Log-Normal, 
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Fig. 10. Average rate of a single-tier cellular network over Log-Normal fad- 
ing {Fg = 1, SNR = 10 dB). Markers show Monte Carlo simulations. Solid 
lines show the analytical framework, which is computed by using Theorem\l\ 
Corollary ^v/ith e = 0.05, Proposition^wiXh = - in (10) 0-^/20 dB 
and Nqhq = 5, and Remark [2] with Nqqq = 2000. Furthermore, 
Mo (•) = Mi (■) are obtained from [75, Eq. (2.54)] with iVoHQ = 5. 




Fig. 11. Average rate of a single-tier cellular network over composite 
Nakagami-m and Log-Normal fading (Fg = 1, SNR = 10 dB, and 
A = 0.25). Markers show Monte Carlo simulations. Solid lines show the 
analytical framework, which is computed by using Theorem [T] Corollary [2] 
with £ = 0.05, Proposition^with fi = - \n (10) o-^/20 dB and Nqhq = 5, 
and Remark^vjith Nqcq = 2000. Furthermore, Aio ( ) = -Mi ( ) are 
obtained from [75, Eq. (2.58)] with Nqhq = 5. 



and composite Rice and Log-Normal distributions, as well 
as that it increases for less severe fading (m increases) and 
in the presence of a stronger line-of-sight component (K 
increases); and ii) the average rate strongly depends on the 
shadowing standard deviation a of Log-Normal, composite 
Nakagami-m and Log-Normal, and composite Rice and Log- 
Normal distributions, as well as that it decreases significantly 
for more severe shadowing (a increases). 

d) Correlated Log-Normal Shadowing: In Figs. [T3] and 
[l4l numerical examples in the presence of shadowing corre- 
lation over a composite Nakagami-m and Log-Normal fading 
channel are shown. We observe that the proposed approach 



Fig. 12. Average rate of a single-tier cellular network over composite Rice 
and Log-Normal fading (Fg = 1, SNR = 10 dB, and A = 0.25). Markers 
show Monte Carlo simulations. Solid lines show the analytical framework, 
which is computed by using Theorem [T] Corollary |2] with e = 0.05, 
Propo.«Y!on |4] with p, = — In (10) (T^/20dB and A^ghq = 5. and Remark 
12] with Nqcq = 2000. In Proposition O the series in I17t is truncated 
to the first 100 terms. The application of Remark [3] provides the same 
result and accuracy, but with less computational complexity. Furthermore, 
Mo (■) = Mi (■) are obtained from l l46t with A^ghq = 5. 
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Fig. 13. Average rate of a single-tier cellular network over correlated 
composite Nakagami-m and Log-Normal fading {Fb = 1). Markers show 
Monte Carlo simula tion s. Solid lines show the analytical framework, which is 
computed by using )221 with A^ghq = 5, Theorem\l\ CoroZ/ary [2] with e = 
0.05, Proposition^with {m = 2.5, a = 6dB, p. = - In (10) 0-2/20 dB, 
a = 4) and A^ghq = 5, and Remark^v/ith Nqqq = 1000. Furthermore, 
Mo (■) = Mi (■) are obtained from [75, Eq. (2.58)] with A^ghq = 5. 



for correlated Log-Normal shadowing is very accurate for 
different choices of the correlation coefficient. Furthermore, 
by comparing Figs. [13] and [14] with Fig. IH we note that the 
framework in Section IIII-FI reduces to the independent case 
for p — 0. More specifically, the figures show the following 
performance trends: i) for high-SNR, the average rate is 
independent of the density of BSs regardless of the shad- 
owing correlation coefficient; ii) for low-SNR, the average 
rate slightly decreases if the shadowing correlation coefficient 
increases; and iii) for high-SNR, the average rate increases if 




Fig. 14. Average rate of a single-tier cellular network over correlated 
composite Nakagami-m and Log-Normal fading {Fg = 1). Markers show 
Monte Carlo simula tion s. Solid lines show the analytical framework, which is 
computed by using )22t with Nqhq = 5, Theorem^ Coro/Zan' [2] with e = 
0.05, f ropoiiYion [3] with (m = 2.5, cr = 6 dB, = - In (10) O-2/20 dB, 
Cf = 5) and A^ghq = 5, and Remark^with Nqcq = 1000. Furthermore, 
Mo (■) = Mi (■) are obtained from [75, Eq. (2.58)] with AfcHQ = 5. 



Fig. 15. Average rate of a single-tier cellular network with frequency reuse 
over composite Nakagami-Hi and Log-Normal fading (A = 0.25). Markers 
show Monte Carlo simulations. Solid lines show the analytical framework, 
which is computed by using Theorem [T] Corollary |2] with e = 0.05, 
Proposition ^m\h {m = 2.5, cr = 6dB, = - In (10) 0-2/20 dB) 
and Nqhq = 5, and Remark [2] with Nqqq = 2000. Furthermore, 
Mo (■) = Mi (•) are obtained from [75, Eq. (2.58)] with A^ghq = 5. 



the shadowing correlation coefficient increases. Even though it 
may seem counterintuitive that the average rate increases with 
shadowing correlation, this result seems to originate from the 
BSs association policy adopted in the present paper (the probe 
mobile terminal is associated with the closest BS). Finally, 
we emphasize that (apparently) counterintuitive trends in the 
presence of Log-Normal shadowing have been observed in 
other papers, e.g., [22] and [97] where it is shown that the cov- 
erage probability increases and that the blocking probability is 
not always increasing with the shadowing standard deviation, 
respectively. This confirms, once again, the importance of 
taking into account correlated Log-Normal shadowing for 
accurate performance prediction. 

e) Impact of Frequency Reuse: In Fig.fTSl the impact of 
frequency reuse on the average rate of a single-tier cellular 
network over composite Nakagami-m and Log-Normal fading 
is studied. Our numerical analysis confirms the findings in 
Section IIII-EI and that frequency reuse is detrimental for the 
average rate. As suggested in [22], frequency reuse signif- 
icantly improves the coverage probability. The analysis of 
this trade-off is out of the scope of the present paper but 
is currently being investigated by the authors. 

f) Framework Validation for Multi-Tier Cellular Net- 
works: In Fig. [16] and Fig. [17] the average rate of a two- 
tier cellular network over composite Nakagami-m and Log- 
Normal and composite Rice and Log-Normal fading is ana- 
lyzed, respectively. As an illustrative example, we consider the 
situation where the BSs of every tier transmit with a power 
that is inversely proportional to their spatial density. This is 
a reasonable choice if, e.g., the first tier is used to model 
macro BSs and the second tier is used to model femto BSs. 
We observe, as expected, that the average rate significantly 
increases when the BSs density of the lower tier increases. 
The framework provides a very good accuracy. Furthermore, 
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Fig. 16. Average rate of a two-tier cellular network over composite 



(t) 



1 for t 



1,2, and 



Nakagami-m and Log-Normal fading (F^ 
SNR = SNRi). Markers show Monte Carlo simulations. Solid lines show 
the analytical framework, which is computed by using Corollary^ Corollary 
121 with e = 0.05, Proposition [3] with (m = 2.5, a = 6dB, p = 
-ln(10)o-2/20dB) and A^ghq = 5, and Remar/t |2] with Nqcq = 2000. 
Furthermore, Mt o (■) = Mt i (■) for t = 1,2 are obtained from [75, Eq. 
(2.58)] with Nqhq = 5. 



we note a negligible difference between the two fading models 
for the chosen set of parameters. 

g) High-SNR Scenario: Finally, we observe that, in all 
figures, the average rate increases with the SNR tending 
towards a horizontal asymptote for high-SNR. By direct 
inspection, the reader can verify that the horizontal asymptote 
coincides with the average rate that can be computed by 

' (A„„) 



(SNR, 



= 7^ 



using Corollary [3] and Corollary [8] with TZ] 
in Corollary [9] In other words, the SNR region where the 
average rate is flat corresponds to the interference-limited 
operating regime. As a consequence, we conclude that the 
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Fig. 17. Average rate of a two-tier cellular network over composite Rice and 
Log-Normal fading (F^*' = 1 for t = 1,2, and SNR = SNRi). Markers 
show Monte Carlo simulations. Solid lines show the analytical framework, 
which is computed by using Corollary |6] Corollary |2] with e = 0.05, 
Proposition^viWh (K = 10, cr = 6 dB, = - In (10) 0-2/20 dB) and 
A^GHQ = 5, and Remark |2] with Nqcq = 2000. In Proposition |3] the 
series in 1171 is truncated to the first 100 terms. The application of Remark\3\ 
provides the same result and accuracy, but with less computational complexity. 
Furthermore, Mt,o ( ) = -Mtj (■) for i = 1, 2 are obtained from l l46t with 
A^GHQ = 5. ' 



simple formulas given in these corollaries are quite accurate 
for SNR > SNR*, where SNR* is the corner point where the 
average rate starts approaching the horizontal asymptote. For 
example, SNR* w 20dB for a = 4 in Fig.[I] In general, SNR* 
depends on the density of BSs, the path-loss exponent, and 
the fading channel model. The reader can identify SNR* for 
different cellular setups by direct inspection of all the figures 
shown in the paper Figure |8] shows many cellular networks 
deployments of interest. 

By observing the asymptotic behavior of the average rate 
for high-SNR, an interesting problem to investigate is whether 
SNR* is lower or greater than the typical operating SNR of 
current cellular networks deployments. This question is inter- 
esting because many papers assume, for analytical tractability, 
that the additive noise is always negligible and that, as a 
consequence, cellular networks are interference-limited. All 
the figures shown in our paper confirm that the interference- 
limited assumption is accurate only if the operating SNR is 
greater than SNR*. In order to assess whether the operating 
condition SNR > SNR* is usually verified, we investigate, 
as an example, the same setup as in [43], which holds for 
typical cellular networks. More specifically, we consider: i) 
SNR = P/<7%; ii) a% = W/Lq; iii) W = fcsToBo = 
— 104dBm — 4 • 10~^'*Watt is the noise power, where 
ks = 1-38 • 10~^^ Joule/Kelvin is the Boltzmann's constant. 
To = 290Kclvin is the noise temperature, Bq = lOMHz is the 
receiver bandwidth; and iv) Lq — — 38.5dB = 1.41 • 10^^ is 
the path-loss at a reference distance of one meter Accordingly, 
we have SNR = P/a% = {Lo/W)P = 3.52 • lO^P. By 
assuming typical transmit powers equal to PMacro — 40Watt, 
-PMicro = 6.3Watt, Ppico = 0.13Watt, and Ppcmto = 
O.OSWatt for macro, micro, pico, and femto BSs [42], [101], 



we obtain SNRMacro = 111.50dB, SNRmicio = 103.45dB, 
SNRpico = 86.60dB, and SNRpcmto = 82.45dB, respec- 
tively. By comparing these typical operating SNRs with the 
SNR* shown in our figures, we conclude that the condition 
SNR > SNR* is well satisfied for many system setups 
analyzed in the paper For example, let us consider the setups 
shown in Fig. E] We observe that SNR > SNR* is always 
verified in dense (A = 10^^) cellular networks for every path- 
loss exponents, and in medium/sparse (A = 10^'*) cellular 
networks for low path-loss exponents (a < 4 in the figure). 
On the other hand, for larger path-loss exponents the condition 
may not be verified for some types of BSs. This confirms that 
typical cellular networks deployments can be approximated, 
in most cases, to be interference-limited, and that the simple 
frameworks in Corollary |3] and Corollary [8] can be used for 
first-order performance analysis, design, and optimization. 



VI. Conclusion 

In this paper, we have introduced a comprehensive math- 
ematical framework for the analysis of the average rate of 
multi-tier cellular networks whose BSs are assumed to be 
randomly distributed according to a PPP spatial distribution. 
The framework is applicable to general fading channel models 
with arbitrary fading parameters. Numerically efficient and 
stable algorithms to compute some transcendental functions, 
such as the Meijer G-function, have been proposed. The 
framework needs the computation of either single- or two-fold 
integrals for general fading distributions and arbitrary path- 
loss exponents. Furthermore, shadowing correlation can be 
taken into account with another extra numerical integral. The 
framework can handle random frequency reuse, and it simpli- 
fies significantly for interference-limited cellular networks and 
for high-SNR setups. Extensive Monte Carlo simulations have 
confirmed the accuracy of the proposed analytical methodol- 
ogy for various fading distributions and cellular deployments. 



Appendix I 
Proof of Theorem\l\ 

By using [66, Lemma 1] with N = M = 1, we have: 

X 

''+°° My (z) - Mx,Y {z) 



E On 1 



■exp{— z}(iz (32) 



(a) 



+00 



My{z) [l~Mx{z)] 



exp {—z} dz 



where: i) X and Y are arbitrary non-negative random vari- 
ables; ii) Mx,Y (z) = E {e" 



^z{X+Y) 



} is the MGF of random 
variable X + Y; and iii) (a) holds if X and Y are independent. 

From m with X = (P.goe"")/<, Y = I^g^iO/'^ 



'N' 



and by taking into account that, conditioning upon X and 
Y are independent, the expectation in (I?) can be re-written as 
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n A) ^ ^ / ^ exp {-.AZ. (.) C} exp {- } 

■ C^-i exp {-^\Zj (z) a (I) [ttAZ, (z)]-^ 



(a) 1 
< 



(38) 



SNR, 



follows: 



E <^ In 1 



< + /agg (0 

expj— 2| 



Mi^^^ (z; [1 - Mo (SNRr"^)] dz 

(33) 

We note that the identity in ( |33] | avoids the need of 
computing Pcov and makes our analytical development totally 
different from current practice [20], [22], and [43]. 

From (|33] |. it is apparent that a closed-form expression 
of Aii^ {')£,) is needed. This is the MGF of the aggregate 
interference, which is generated by all the interferers that lie 
outside a disk of radius In other words, due to the tier and 
BS association policy, MTq has an exclusion zone around 
it where no interfering BSs are located and are allowed to 
transmit. This is called exclusion region [5]. The MGF of 
^agg (•) has been studied in [91] for a generic annular region 
with radii A and B. In particular, Mj^ (•; can be obtained 
from [91, Eq. (6)] by letting A ^ ^ and B — > +c», as follows: 

X exp|-^As2/"r (^1 - E{g,'/"} 

X exp I^ttXs^^E |5,'/"r (^1 - ^, sfftr" 

(34) 

In [91], the expectation over the fading distribution of the 
interference channels is not computed in closed-form and only 
bounds for Rayleigh fading are provided. To the best of our 
knowledge, there is no analytical framework that provides an 
exact and closed-form expression of the MGF in ( |34] i for 
general fading channels. In what follows, we provide a general 
methodology to this end. This is a contribution of this paper. 

By using [69, Eq. (6.5.3)], [69, Eq. (6.5.4)], and [69, Eq. 
(6.5.29)], we have: 



r (z, x) — r (z) — 7 (z, x) 

= r (z) - r (z) x^ exp {-x} 



+ 00 



r(z 



1 



(35) 



(36) 



By substituting in -^/agg (sf) simplifies to: 

-^/.gg (s; = exp {nXe } exp {-nXeMi (^r") } 

xexp{-TrXeTi{sr'')} 

where 7/ (•) is given in (|9|l. Closed-form expressions of 7} (•) 
are available in Propositions\l\^ 

Finally, by substituting ( [36l ) in ( [33] ), the average rate in O 
simplifies to (O by using the change of variable y = z^^" 



and by applying the integration by parts to the integral in ^. 
This concludes the proof. □ 



Appendix II 
Proof of Remark\l} 

From (ISJ, by using the change of variable SNRj/ — z, we 
have: 

+°°l-Mo{z)dz a f+°° 1 - Mo (z) , 



n= \' ' — -/ ' n{z-x)dz 

Jo \z) z 2 Ja Zi z 

(37) 

where TZ{-;-) is given in ( |38] | at the top of this 
page, and (a) follows by taking into account that < 
exp {- (z/SNR)^"/2} < 1 for every choice of the param- 
eters, while (b) follows from [73, Eq. (2.1.1.1)]. 

As a consequence, the average rate in (|J7] | is upper- and 
lower-bound as follows: 

l^Moiz}dz_- 1 - Mo (z) dz 

^ -Jo Zi{z) z 

(39) 



10 Zi {z) z 

where: 



(f)<'^)- 



SNR 2 

[l-XoW]2/^+'^ (z) dz 



(40) 







Equation ( fTOb immediately follows from (|40] | by observing 
that limA^+oo (A) = for every a> 2. This concludes the 
proof. □ 



Appendix III 
Proof of Proposition^ 

By using [75, Eq. (2.21)] and [73, Eq. (2.2.1.2)], Mf^ (•) 
in (|9]l can be computed as follows: 

(41) 

By substituting W\\ in 7} (•) in (|9l) and by using the identity 
[69, Eq. (15.1.1)]: 



r (fc + A) 5^+*^ 

■^^T{k + B){s + Cr'-'' 



(42) 



C 



{s + cyT{By 

with A, B, C, X, and s being positive constants, we eventually 
obtain ( fT4] i with some algebraic manipulations and using the 
identity F (z + 1) = zF (z). This concludes the proof. □ 
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fg, (x) « (1 + exp {-K} — = WnCjn exp {- (1 + K) Cb^x} h (2^ K {I + K) Cor,x^ 

^ n—1 



V n=l 



+ K+ [s/Ojr,) 



exp < — s 



1 + K+ (s/w„) 



(46) 
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(47) 
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^^GHQ 

E 



+00 



i)„a)„ y 



r (2 + fc) 



- l2+fe- 



1^1 fc + 2,1, 



if(l + if) 

1 + if + (s/Wn) 



(48) 



Appendix IV 
Proof of Proposition^ 

By using the approximate expression of the PDF of Log- 
Normal random variables in [77, Table IV], i.e., fg^ (x) « 
E^fr ^nS [x - io(^-"+A')/io^, Xf) (.) in @ 
can be computed as follows: 



'10 



n=l 



(43) 



X exp 



|_-^q(x/2^S„+m)/10^| 



By substituting ( |43] | in 7/ (•) in (|9]l and by using the identity 
[69, Eq. (13.1.2)]: 



+°° gX + k 



E 

fc=0 



T{k + B) r{B) 



(44) 



with B and z being positive constants, we obtain (fTsl l after 
some algebra. This concludes the proof. □ 

Appendix V 
Proof of Proposition^ 

Similar to Appendix |IV] in order to have an analyti- 
cally tractable expression of PDF and MGF of composite 
Nakagami-m and Log-Normal fading, we use the approxi- 
mation of the Log-Normal distribution that is based on the 
Gauss-Hermite quadrature in [75, Eq. (2.58)] and [76, Table 
V]. Accordingly, A^f ^ (•) in © can be computed as follows: 

M} ' (s) 



/tt r (m) 

-'Vghq 



f 1) 
(45) 



By substituting ( l45T l in 7/ (•) in ©> from (l42l i we eventually 
get ( fTSI l. This concludes the proof. □ 

Appendix VI 
Proof of Proposition\4\ 

Similar to Appendix |iy] in order to have an analytically 
tractable expression of PDF and MGF of composite Rice 
and Log-Normal fading, we use the approximation of the 
Log-Normal distribution that is based on the Gauss-Hermite 
quadrature in [77, Table IV]. In particular, from [78, Eq. (6)] 
we obtain ( |46] ) shown at the top of this page, where (a) is 
obtained by using [75, Eq. (2.17)]. Thus, A^^^ (•) in ^ can 
be computed as shown in ( |47| | at the top of this page, where 
(a) is obtained by using [69, Eq. (6.631)]. By substituting ( |47| | 
in 7} (•) in we obtain (l48l l shown at the top of this page. 

Since the infinite series in (|48] | is not fast converging, we 
elaborate further 7} (•) in order to obtain a better expression 
for simple numerical computation. For if 7^ 0, this can be 
obtained by first replacing iFi {■,-,■) with its series expansion 
in [69, Eq. (13.1.2)], and then computing the infinite series in k 
by using (l42l l. Eventually, we obtain ( fTTl ) with some algebraic 
manipulations and using the identity T {z + 1) = zT (z). The 
case if = can be obtained by noting that iFi {A, 1,0) — 1 
for every A, and by applying the same procedure as in 
Appendix FV] with m = 1. This concludes the proof. □ 

Appendix VII 
Proof of Corollary^ 
Let z^(asymptoto) ^ lim^^o+ W (z) in (HHl. This Hmit 
can be computed as follows: 

^(asymptote) _ ,;™ /^Qd,Qn („ ^(o^-ZVi 



(a) ,. 
~ lim 

C-i-+oo 



1 



A(ai3,0) 
1- A(ai,,0) 



(49) 
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where C = 1/z and (a) is obtained by using [68, Eq. 
(8.2.2.14)]. 

From (gill, (HHl can be obtained by using [92, Theorem 
1.8.3], and more explicitly by using [92, Eq. 1.8.8] whose 
parameters are defined in [92, Eq. 1.4.2] and [92, Eq. 1.4.7], as 
well as by taking into account that lim2_^Q+ pFq (a, b, z) — 1, 
and, thus, the generalized hypergeometric function pFq {■,-,■) 
can be neglected. This concludes the proof. □ 
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